

%% load estimation results
load('cold_war_replication.mat')
[~,use] = min(l0);
thetaStar = theta(:,use);
clearvars -except thetaStar K mpow N nmissing reb X Y YY Z specif


%% load everything of interest
load('../../baseline_model/sim_eq_round2.mat')
N=sum(nmissing);
intervEq=intervEq(:,nmissing);
U0=U0(:,:,nmissing);
V_R=V_R(:,:,nmissing);
V_US=V_US(:,:,nmissing);
V_UK=V_UK(:,:,nmissing);
V_France=V_France(:,:,nmissing);
V_Russia=V_Russia(:,:,nmissing);
V_China=V_China(:,:,nmissing);
[DMES,intervEq] = createDMES(intervEq,V_R,V_US,V_UK,V_France,V_Russia,V_China,YY,specif);


%% standard errors via Hessian
InvSigmaH = HessPhat(thetaStar,1,Y,X,Z,YY,V_R,V_US,V_UK,V_France,V_Russia,V_China,U0,intervEq,DMES,reb,mpow);
seH = sqrt(diag(inv(InvSigmaH)));
resultsH = printResults([thetaStar,seH,2*normcdf(abs(thetaStar./seH),'upper')],reb,mpow,specif);
disp(' ')
disp('for Tables E1-E2:')
disp(resultsH)


%%
clearvars -except thetaStar InvSigmaH seH resultsH
save('inference_cold_war.mat')



